Thermal Rectification in Billiard-like Systems 
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We study the thermal rectification phenomenon in "billiard" systems with interacting particles. This interac- 
tion induces a local dynamical response of the billiard to an external thermodynamic gradient. To explain this 
dynamical effect we study the steady state of an asymmetric billiard in terms of the particle and energy reflec- 
tion coefficients. This allows us to obtain expressions for the region in parameter space where large thermal 
rectifications are expected. Our results are confirmed by extensive numerical simulations. 
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Introduction. The emergence of macroscopic behavior from 
the underlying microscopic dynamical laws is a major puz- 
zle of theoretical physics. However, despite its long history, 
attempts to link these macroscopic laws to the underlying mi- 
croscopic dynamics have not been conclusive thus far [1]. 

With the advent of nanotechnology, the study of dynamics 
of systems at a nanoscopic scale starts to play a central role. 
In view of the broad spectrum of possible applications, the in- 
vestigation on how the properties of the microscopic dynamics 
can be used to understand and thus, to control the behavior of 
nano-devices has attracted much attention [2]. 

One such example is the phenomenon of thermal rectifica- 
tion (TR). A thermal rectifier is a device in which the magni- 
tude of the heat current depends on the sign of the imposed 
temperature gradient. In the past there have been several ex- 
perimental attempts to design such a device (see e.g., [10] and 
references therein). However, a possible theoretical mecha- 
nism to explain the phenomenon from first principles has been 
proposed only recently [6]. 

So far, there is no unique framework to explain the phe- 
nomenon of thermal rectification. In [6] and the later papers 
[7, 8], different models of anharmonic oscillators coupled to 
stochastic heat baths have been studied. Because of anhar- 
monic terms, the width and position of an "effective" phonon 
band depend on temperature. It has been shown that by tuning 
the anharmonicities of two or more different coupled chains it 
is possible to control the effective overlap between the phonon 
bands of the two chains. In particular, the extent of this over- 
lap can be made to depend on the sign of the imposed thermal 
gradient, thus leading to a TR. A further step to devise a ther- 
mal diode has been discussed in [8] in terms of the negative 
differential thermal resistance observed in some anharmonic 
chains. More recently, in [9] a thermal rectifier based on a 
spin-boson nanojunction model has been discussed. 

In this Letter we study the phenomenon of thermal recti- 
fication in low dimensional "billiard" Hamiltonian systems. 
Roughly speaking, a billiard consists of an ensemble of free 
particles colliding elastically with a fixed boundary. Further- 
more, to study heat transport and thus rectification, we con- 
sider that the billiard is coupled at two different openings of 
its boundary with two different thermodynamic baths. 

In a billiard, as described above, TR is identically zero irre- 



spective of the geometry of the boundary, as the strengths of 
the fluxes inside the billiard are determined only by its cou- 
pling with the external baths. Here, we show that if the parti- 
cles are allowed to interact (as defined below), then rectifica- 
tion of the heat flux is possible, and can be quite large. 

The model we consider consists of a gas of non- 
interacting point particles of mass m that move freely in- 
side a one-dimensional open channel made of 2 identical two- 
dimensional cells (like chaotic billiards), connected through 
openings of size S. At the center of each cell, there is a fixed, 
freely rotating disc of radius R with which particles interact, 
exchanging energy as follows: when a particle hits a disc, it 
exchanges its tangential velocity v t with the angular velocity 
oj of the disc, while the normal component v n is reflected [11]: 



(1) 



This model was introduced in [11, 12] (with a different ge- 
ometry) to study macroscopic transport in terms of its micro- 
scopic dynamics. It was found that the disc mediated an ef- 
fective interaction among the particles that leads to a realistic 
macroscopic transport behavior. Moreover, at each cell this 
effective interaction is sufficient for the establishment of ap- 
proximate local thermal equilibrium [11, 13]. 

As a result of the interaction, when a particle enters the bil- 
liard from one of the reservoirs, its probability to reach the 
other reservoir, i.e., the transmission coefficient, depends not 
only on the geometry of the billiard but also on the local ther- 
modynamic gradients. This dynamical effect was studied in 
detail in [14] in a billiard of (effectively) interacting particles 
with a geometry similar to the one that we describe below (see 
also [13]). 

To set up the thermal rectifier, we exploit the dependence of 
the dynamical effect described above on the geometry of the 
cell. We choose the size of the disc in the left cell (-) smaller 
than the disc in the right cell (+). The transmission coefficient 
will then depend differently on the thermodynamic fields, and 
this leads to TR. The geometry of our setup is described in 
Fig. 1. 

To impose nonequilibrium, the two cells are connected at 
the ends to two reservoirs of particles through openings of the 
same size 5. The reservoirs are idealized as infinite chambers 
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FIG. 1 : Geometry of the rectifier. The boundary of each cell is made 
of four arcs of circles of radius R s . The radius of the rotating discs is 
chosen to be larger than the width of the openings 5, so that no par- 
ticle can cross the cell without undergoing any collision. The arrows 
represent the flows of particles across the openings and the symbols 
fix the notation for the rates of particles j and of energy q. Unless 
otherwise specified, the numerical results in this letter correspond to 
a rectifier with discs of radius 0.04 for the (-) cell and 0.18 for the 
(+) cell, 5 = 0.04 and R s = 1. 



containing an ideal gas at a certain density n and temperature 
T[12]. 

From the left reservoir, particles are injected into the (-) cell 
at a constant rate j > with Boltzmann distributed energy with a 
mean temperature T > . Analogously, the right reservoir injects 
particles into the (+) cell at rate j < and temperature T < . The 
rate at which energy is injected into the system from the left 
and right reservoirs is then given by q > = *j > T > and q < = 
^j < T < respectively [13]. Each injected particle is eventually 
re-absorbed by one of the reservoirs and this happens when it 
hits one of the reservoirs. 

The density n and temperature T of the reservoir are related 
with its injection rates as 

'3/2 
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At the junction between the two cells particles cross from 
(-) to (+) at a rate j R and from (+) to (-) at a rate j L . We denote 
by q R and q L the corresponding energy rates. 

In the steady state the conservation of particle number and 
of energy can be expressed as a set of balance equations, [14] 



Jl = 



{ 1 - a h)j> + ( 1 - a R)j< 
1 - a l a R 

7^(1-7l)<7> + (1-7r>< 



(2) 
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where is the reflection probability for the particles that are 
injected into the cell (-) from the left and a R for the particles 
injected into the cell (-) from the right and similarly for the (+) 
cell. The parameters 7 stand for the energy reflection coeffi- 
cients with the same indexing as the a's. Analogous equations 
are obtained for j R and q R , by exchanging L <-> R, + <-> — , 
and > <-> <. 



Note that out of equilibrium the reflection from the left is 
not necessarily the same to the reflection from the right due to 
the effective particle interaction. 

In terms of the rates at the junction the particle current is 
obtained by $ = j R — j L and analogously, the heat current 
is <£> M = q R — q L . Using (2) the currents are 
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(3) 



where = (l — xZ ) (l—Xy), with x — {a, j},y = 
{L, R}, are the total transmission probabilities for the par- 
ticles that enter the channel from the y side. The currents 
in (3) are obtained as the difference of the products of these 
transmission probabilities with the injection rates, weighted 
by combinatorial factors that account for multiple reflections 
between (-) and (+). 

Reflection coefficients for a single cell. In [14] a phe- 
nomenological theory for the reflection probabilities <x^ was 

R 

obtained. These results are for a single cell with injection of 
particles and energy at rates j > , g > from the left and j <t q K 
from the right. The reflection probabilities depend on the local 
thermodynamic fields, that is 



\3>,Q.»3<,<1<} 57^ =7 fe {3 > ,q > ,J < ,q <J 



(4) 



Since the dynamics is homogeneous in energy, the number of 
independent variables in (4) is reduced to three. It was found 
convenient to distinguish two contributions: one (called a G 
and 7 G ), of purely geometrical origin and another (called e 
and v), of dynamical origin arising as the system is driven out 
of equilibrium: 

«l =a G (n) + £l (n, Sj, 5q) , 
7l =7g ( n ) + ^ (n,6j,6q) , 

R R 

where Sj — (j> — J<)/(j> + 3 <) an d similarly for 8q. By 
definition, e must vanish at equilibrium, i.e., e L (n, 0, 0) = 0, 

R 

and by left-right symmetry of the single cell, e L = — e R . 
Detailed analysis showed that: 



At equilibrium a L 



,(n). 



• a G (n) saturates, for low and for high densities, to two 
constant values connected by a linear dependence for 
intermediate densities. 

• £ vanishes on a submanifold of (J > , g > , j < , q K ) where 
T > = T < . This means that deviations from the equilib- 
rium a 's only occur for Sj 0. 

• Close to equilibrium, the dynamic component depends 
linearly on the local gradients of the particle injection, 
namely e L = ±A5j, = ±B~Sj Note that v depends 

R. R 

on Sj, not on 8q. 
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FIG. 2: a) Geometrical reflections ao(n) and 7g(^) for a single cell 
with the same dimensions as cell (+), maintained at equilibrium, b) 
Dynamical reflections £l and i^l as a function of the local gradients 
of particle injection rates for the cell (+). Note that by definition 
ttG = (ql + Qr)/2 and similar for 7g- In both panels, the symbols 
correspond to the energy reflection (hashed circles) and to the particle 
reflection (open squares) and the dashed lines correspond to linear 
fits. 



In Fig. 2, we show the behavior of these two contributions 
for a cell of type (+). These results were obtained under very 
general reasoning and thus, they are expected to hold quali- 
tatively for different models of interactions as long as the dy- 
namics inside each cell is sufficiently mixed and phase-space 
volume is preserved [14]. 

We have verified that this qualitative behavior is the same 
for different sizes of the disc, the only change being the val- 
ues for the saturation plateaus of the geometric contribution 
(Fig. 2-a) and the slope of the linear dependence of the dynam- 
ical contribution on the local gradients (Fig. 2-b). In particular 
we have found (not shown) that the variation, from the lower 
to the upper plateau, of a G is much smaller (~ 20 times) for 
the (-) cell than for the (+) cell and similar for 7 G . Therefore, 
these quantitative values characterize the geometry of the cell 
and its dynamical response to the presence of thermodynamic 
gradients. 

Combining the two cells. We now apply these formulas for 
the two-cell case: For the (-) cell, the particle currents are j > 
and j L , while for the (+), they are j R and j < . From our find- 
ings for a single cell, we find, to a very good approximation: 



R J>+Jl 



> C L 



_+ _ ^+ 3r 3< 

" R JR+J< 



(5) 



3> + Jl 



3 R +3< 



(6) 

Here and below ± refers to one of the two cells. Substituting 
(5) into (2) we obtain the currents j L at the center of the cell in 

R 

terms of the parameters that define the geometry of the cells 
as 



Jl = 
3r = 



Qti 1 ~ Qa)j> + (! - et)j< 
1 - QaBt 

(l-gq)j> +gg(l~gq)i< 

1 - QaQt 



(7) 



with effective reflection coefficients given by — a G (n )— 
A and with n the mean density in the respective cell. One 
finds finally for the particle current (to within our approxima- 
tion) 



(1 - QZ) (1 - 9 + a ) 
1 - Q a Qt 



\3< - .7; 



(8) 



Thus, near equilibrium, <& n does not depend on the energy 
field q. In contrast, from (6) it is clear that the current of 
energy depends on both fields. Substituting (7) into (6) we 
find for the dynamical reflections 



S" 



B 



2(l-g-g+)j> 

(l-«+)(i<-i>) 



2(i-g-e+))j < 
(i-«;)(j'<-i>) 



(9) 

Finally, inserting (6) and (9) in (3), one obtains an expression 
for the energy current exclusively in terms of the parameters 

®u = ®u ti> , 9> . 3< . 9< ; et , 7g > b± )■ 

Inverting the imposed fields. From now on we refer to the 

case in which the injections at the left are (j > , g > ) and at 
the right (j < ,q K ) as positive gradient (PG) and to the in- 
verted case as negative gradient (NG). Due to the effective 
interaction among particles the local fields in the cells for 
the NG are different from those of PG. Intuitively, this is be- 
cause the temperature at the junction of the two cells changes 
when one goes from PG to NG. In particular, the energy re- 
flection coefficients in (9) are not invariant with respect to 
the sign of the gradient. Let n ± be the particle density for 
PG and let n ± be the particle densities for the NG. Defining 
Qa = q g(" ± ) ~ anc l usm g me tilde quantities in (9), one 
obtains in analogy with (3) for the NG: 



~7 

t l<1< 



-1 



1 



7 L 7 R 



7± ± P± 



with rl = (1 - 7 L )(1 - R ) and 7^ 

R R R 

To quantify the TR we define the thermal rectification index 
as the quotient between the heat currents for the PG and NG 



e„ = 



max{|$J,|$J} 
min{|$J, |*J} 



so that G n > 1. For zero rectification <d u = 1, 
Results. In Fig. 3, we show in a density plot the logarithm of 
the thermal rectification index 9 M as a function of the external 
gradient. A large TR is found in a region where Sj and Sq are 
linearly related. The maximal value of 9 U observed in Fig. 3 
corresponds to heat currents $ u and $ u where one is three 
orders of magnitude smaller than the other. 

As found in previous models, the magnitude of TR depends 
on the external gradients. For our type of models the largest 
values of G u are obtained when either $ u or <I> U (but no both) 
is very small. This means that TR is maximal when for a 
given external gradient the system behaves as an insulator. If 
for a given gradient (6j, Sq), $ u = then for the inverted 
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FIG. 3: Color density plot of the logarithm of the thermal rectifica- 
tion index log (Q u ) as a function of the external gradients Sj and 
5q for the billiard described in Fig. 1. By construction the data is 
symmetric along the diagonal. The dashed line correspond fit to the 
numerical data of the zero-current condition (10). 



gradient (Sj, —Sq), $ u will not in general be zero. This is 
because the reflection probabilities will be different through 
their dependence on the local thermodynamic fields. 

Equating <f> M to zero we obtain a condition that involves the 
reflection probabilities 7: 

(1-7l)(1-7l-) = flm 
(1-7r)(1-7£) 1> ' 

This equation states the values of the external gradient for 
which our model behaves as a heat insulator. The dashed line 
in Fig. 3 corresponds to a linear fit to the condition (10) com- 
puted from the values obtained from the numerical simulation. 
Discussion and conclusions. In our model TR is due to the 
dynamical response of each cell to a given external gradient 
and thus, to the interacting character of the particles. If both 
cells and their disks are identical, then the response for the PG 
is identical to the one for the NG and there is no TR. On the 
other hand, it is tempting to consider the opposite case when 
the domains are equal, but the cell (-) has no disc at all. While 
we have observed TR for this case, we do not report further 
on this as our analytical treatment is not longer valid. 

As we have seen before, rectification has a geometric com- 
ponent (occurring when Sj = 0). For the billiard of Fig. 1 a 
thermal rectification of A u ks 1.04 is obtained when the geo- 
metric component for the PG and the NG takes values in the 
lower and upper plateaus respectively. 

This rectification power can be enhanced to much larger 
values when the billiard behaves as an insulator for positive or 
negative gradient, and this can be achieved when Sj 7^ 0. 



The geometric component also allows a rectification of the 
particle current. However, this rectification cannot be en- 
hanced as from (8) it is clear that the particle current $ n is 
zero only when Sj = 0. 

For open billiards like the one studied here, the finite parti- 
cle current contributes to the TR only through the dynamical 
component . In practice, this situation corresponds to the 
transport in porous materials where the energy is carried by 
the particles. However, for the analogous closed billiard for 
which only energy is exchange between the system and the 
reservoirs large thermal rectification have been numerically 
observed. 

In this Letter, we have established the possibility to observe 
thermal rectification in billiard-like mechanical systems. We 
have shown that when the gas of particles inside the billiard 
interact, the response of the system to an external thermody- 
namic gradient is dynamically controlled. We have studied 
this dynamical effect in terms of the particle and energy re- 
flection coefficients for which an analytical phenomenological 
formulation in terms of the local thermodynamic gradients is 
possible. This dynamical response leads to possibility that the 
system behaves as an insulator for one gradient and as con- 
ductor with the inverted gradient in which case, large thermal 
rectification is obtained. Our results lead to conjecture that 
TR would be observable in nanoscopic billiard-like systems 
for which the gas of electrons is not free. 
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